
#include "SceneManager.h"
float mypow(float x, float y) {
	return (x<0 && abs(y)<1)? -pow(-x,y): pow(x,y);
}

void func() {
	float t1,t2,t3,t4;
	float A = 1;
	float B = -.25;
	float C = -.85;
	float D = 1.45;
	float E = -4.35;
	A = 1;

	float qf = C - (3*B*B/8);
	float qg = D + (mypow(B,3)/8) - (B*C/2);
	float qh = E - (3*mypow(B/4, 4)) + B*B*C/16 - (B*D/4);
	float cc1 = qf/2;	//cubic equation coefficient 1
	float cc2 = (qf*qf-4*qh)/16;
	float cc3 = -qg*qg/64;
	//solve cubic equation
	int numSolution = 0;
	float y1,y2,y3;
	float i2, i3;//imagin part of y2 and y3
	float cf = (3*cc2-cc1*cc1)/3;
	float cg = (2*cc1*cc1*cc1 - 9*cc1*cc2 + 27*cc3)/27;
	float ch = mypow(cg/2,2) + mypow(cf/3, 3);
	if(ch>0) {//only one real root
		float cR = -(cg/2) + sqrt(ch);
		float cS = mypow(cR, (float)(1.0/3.0));
		float cT = -(cg/2) - sqrt(ch);
		float cU = mypow(cT, float(1.0/3.0));
		y1 = cU +  cS - cc1/3;
		y2 = - (cS+cU)/2 - cc1/3;
		i2 = (cS-cU)*sqrt(3.0)/2;
		y3 = - y2;
		i3 = -i2;
		numSolution = 1;
	}
	else if(isEqual(cf,0) && isEqual(cg,0) && isEqual(ch,0)) {//all three roots are equal and real
		y1 = mypow(cc3,float(1.0/3.0)) - 1;
		numSolution = 3;
	}
	else {//three real root
		float ci = sqrt(cg*cg/4 -ch);
		float cj = mypow(ci, (float)(1.0/3.0));
		float ck = acos(-cg/(2*ci));
		float cl = -cj;
		float cm = cos(ck/3);
		float cn = sqrt(3.0)*sin(ck/3);
		float cp = -(cc1/3);
		y1 = 2*cj*cm+cp;
		y2 = cl*(cm+cn)+cp;
		y3 = cl*(cm-cn)+cp;
		numSolution = 3;
	}
	//back to quartic equation
	if(3 == numSolution) {
		float qp, qq;	//let qq and qp be square root of any two non zero from y1, y2 and y3
		if(isEqual(y1,0)) {
			qp = sqrt(y2);
			qq = sqrt(y3);
		}
		else {
			qp = sqrt(y1);
			qq = sqrt((isEqual(y2,0)?y3:y2));
		}
		float qr = -qg/(8*qq*qp);
		float qs = B/(4*A);
		t1 = qq + qp + qr - qs;
		t2 = qq - qp - qr - qs;
		t3 = -qq + qp - qr - qs;
		t4 = -qq - qp + qr - qs;
		printf("%.3f, %3f, %3f, %3f\n",t1,t2,t3,t4);
	}
	else if( 1 == numSolution) {
		float qp, qpi, qq, qqi;
		//qp+Iqpi is square root of imagin root y2
		//qq+Iqqi is square root of imagin root y3
		qpi = sqrt((sqrt(y2*y2+i2*i2)-y2)/2);
		qp = i2/(2*qpi);
		qqi = -qpi;//sqrt((sqrt(y3*y3+i3*i3)-y3)/2);
		qq = qp;//i2/(2*qqi);
		float qr = -qg/(8*(qp*qq-qpi*qqi));
		float qs = B/(4*A);
		t1 = qp + qq + qr - qs;
		t2 = -(qp + qq) + qr - qs;
		printf("%.3f, %3f\n",t1,t2);
	}
}

int main(int argc, char **argv) {
	char **tt;
	tt = new char*[2];
	tt[0] = new char[20];
	tt[1] = new char[20];
	strcpy(tt[0],"assignment4");
	strcpy(tt[1], "spec1.txt");
	//func();
	return SceneManager::my_glmain(argc, tt);
}